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ABSTRACT 


In  this  paper  we  obtain  time-domain  equations  — state  space 
equations  — for  a model  of  a horizontally  stratified  nonabsorbtive 
earth  with  vertically  traveling  plane  compressional  waves.  We 
develop  a complete  model  which  generates  primaries  and  all  multi- 
ples; a primaries  model  which  generates  just  the  primary  reflection 
components  of  a seismogram;  and  a secondaries  model  which 
generates  just  the  secondary  reflection  components  of  a seismogram. 
Additionally,  we  obtain  a complete  (canonical)  decomposition  of  a 
seismogram  into  an  interconnection  of  models  each  one  of  which 
generates  just  the  primaries,  secondaries,  etc.  We  also  show  how 
to  obtain  transfer  functions  from  these  models  by  means  of  a recursive 
procedure.  Finally,  we  show  how  our  normal  incidence  models  can 
be  extended  to  non-normal  incident  waves. 
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I.  INTRODUCTION 
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Ray  theory  models  for  seismic  waves  in  layered  media  may  be 
used  to  generate  synthetic  seismograms,  which  are  useful  for  providing 
rough  data  that  can  be  used  for  such  purposes  as  evaluation  of  new 
deconvolution,  migration,  etc.  techniques.  These  models  are  also  useful 
for  identifying  important  parameters,  such  as  reflection  coefficients.  In 
this  paper  our  attention  is  directed  at  obtaining  time -domain  equations  — 
state  space  equations  — for  a model  of  a horizontally  stratified  nonabsorbtive 
earth  with  vertically  traveling  plane  compressional  waves.  Additionally,  we 
shall  demonstrate  how  to  obtain  transfer  functions  from  the  state  space  model. 

Much  work  has  been  done  in  the  past  (Refs.  1-3,  for  example)  in 
developing  the  basic  equations  for  synthetic  seismograms  for  the  layered 
earth  situation  described  above.  We  shall  discuss  some  of  this  work 
very  briefly  in  Section  HI.  In  order  to  distinguish  between  the  earlier 
work  and  our  work,  we  shall  refer  to  models  from  the  earlier  work  as 
traditional  models,  whereas  we  shall  refer  to  our  models  as  new  models. 

As  a preview,  we  state  the  major  differences  between  the  new  and 
traditional  models:  ( 1 ) the  traditional  models  are  derived  by  examining 

wave  effects  at  the  interface  between  adjacent  layers,  whereas  the  new 
models  are  derived  by  examinir  wave  effects  within  a single  layer;  ( 2. ) 
the  traditional  models  assume  equal  travel  times  in  each  layer,  whereas 
the  new  models  do  not;  and  ( 3 ) the  traditional  models  are  z -transform 
models,  whereas  the  new  models  are  time-domain  state  space  models. 

A system  of  K layered  media  is  depicted  in  Figure  1.  We  adopt  the 
convention  of  calling  the  layer  below  layer  K the  basement.  Each  layer 
is  characterized  by  its  one  way  travel  time,  t\  , velocity,  V^,  and 
normal  incidence  reflection  coefficient  r^  (i  = 1 , 2,  . . . , K).  Additionally, 
interface -0  denotes  the  surface  and  is  characterized  by  reflection  co- 
efficient rQ.  In  Figure  1,  m(t)  and  y(t)  denote  the  input  (e.g.,  seismic 
source  signature  from  dynamite,  aquagun,  etc.)  to  the  layered  earth 


Interface  - 0 ( r^) 


- Interface  - 1 ( r, ) 


\ 


Interface  - 2 ( r^) 
Interface  - 3 ( r^) 


Interface  - (K-2)  (r^ 


Interface  - (K-l) 


Interface -K  (rK) 


Figure  I.  System  of  K layered  media 


system  " which  is  applied  at  interface-O,  and  the  output  (i.  e,  , ideal 
seismogram)  of  the  system  which  is  observed  at  the  surface,  res- 
pectively. 

As  in  Refs.  2 & 3,  we  shall  find  it  convenient  to  draw  ray 
diagrams  with  time  displacement  along  the  horizontal  axis,  so  that 
the  rays  appear  to  be  at  non-normal  incidence  and  so  do  not  overlap 
one  another. 

In  Section  II,  we  derive  state  equations  and  transfer  functions 
for  a complete  model,  complete  in  the  sense  that  it  generates 
primaries  and  all  multiples  for  the  K layered  media.  Connections 
between  the  new  (complete)  model  and  traditional  models  are  made  in 
Section  III.  In  Section  IV  we  show  how  to  modify  the  complete  model 
so  as  to  obtain  a primaries  model  — a state  space  model  which 
generates  only  primary  reflections.  The  extension  of  these  results  to 
a secondaries  model  — a state  space  model  which  generates  only 
secondary  reflections  — is  given  in  Section  V.  A generalization  of  the 
Section  IV  and  V results  is  given  in  Section  VI.  Finally,  we  indicate 
how  to  apply  our  new  models  to  non-normal  incidence  in  Section  VII. 

n.  COMPLETE  MODEL 
A.  Derivation  of  Model 

A state  space  model  for  the  system  of  K layered  media,  depicted 
in  Figure  1,  is  derived  in  this  paragraph  under  the  following  modeling 
assumptions:  (1)  plane  compressional  waves,  (2)  horizontally  stratified 

nonabsorbtive  layers  of  different  thicknesses,  and  (3)  normally  incident 
waves.  The  extension  to  non-normally  incident  waves  is  considered  in 
Section  VII. 


* 


In  a marine  environment,  layer  1 can  be  taken  as  water. 


Figure  2.  Layer-k 


l ; 


To  begin,  we  define  the  states  associated  with  a general  layer, 

the  k**1  layer  (Figure  2),  keeping  in  mind  that  the  logical  choice  of 

states,  in  general,  is  made  to  simplify  the  form  of  a model  and  its 

subsequent  utilization.  For  the  purposes  of  the  present  development, 

it  is  convenient  to  associate  , rk  ^ , and  r^  with  the  kth  layer.  The 

compressional  waves  within  the  k-th  layer  are  identified  by  two  states, 

xU(t)  and  x?(t).  Physically,  these  states  represent  reflections.  State 
k k 

x^(t)  denotes  the  upward  moving  wave  (i.  e.  , ray)  which  has  been 

reflected  off  of  the  bottom  interface  of  layer  k;  whereas,  state  xk(t) 

denotes  the  downward  moving  wave  which  has  been  reflected  off  of  the 

u d 

top  interface  of  layer  k.  To  emphasize  the  fact  that  x^  and  x^  represent 
reflections,  we  depict  this  in  Figure  3,  from  which  we  see  that 


u 

*k 


rk  Uk 


( 1 ) 


and 


d 

Xk  = -r 


k-1  bk 


( 2 ) 


Signal  is  the  result  of  two  phenomena  (Figure  4),  a reflection 
of  an  upgoing  wave  off  of  interface  k-1  — which,  by  definition,  is  x^  — 
and  a transmission  of  a downgoing  wave  through  interface  k-1.  This 
downgoing  wave  is  not  a state  (since  it  is  not  a reflection);  but,  can 
be  expressed  in  terms  of  state  x^  ^ by  the  construction  depicted  in 


-4- 


Figure  4.  Components  of 


^k 


Figure  4.  From  Figure  4,  we  see  that 


In  order  to  complete  our  description  of  x^  , via  Eqs.  (1  ) and  (3)  we 

must  specify  the  functional  dependences  of  x,U  and  u,  on  time.  Observe. 

k k 

from  Figures  3 and  4 that  Eq.  ( 1 ) relates  u and  xu  at  the  bottom  of 

k k 

layer  k,  whereas  w has  been  defined  in  terms  of  x^  and  xU  at  the 

k k k - 1 

top  of  layer  k.  In  order  to  account  for  the  sec.  time  delay  in  traveling 
through  layer  k,  we  write  Eqs.  (1)  and  (3)  together,  as 


-5- 


Figure  5.  Components  of 


u 


,t+Tk» 


rk  *k<t) + rk 


in 


( 4 ) 


Next,  we  direct  our  attention  at  signal  E , which  is  also  the 

k 

result  of  two  phenomena  (Figure  5),  a reflection  of  a downgoing  wave 
off  of  interface  k — which,  by  definition,  is  — and  a transmission 
of  an  upgoing  wave  through  interface  k.  This  upgoing  wave  is  not  a 
state;  but,  can  be  expressed  in  terms  of  x^  by  the  construction 
depicted  in  Figure  5.  From  Figure  5,  we  see  that 


5k  = ^ ' (^)  Xk+1 


( 5 


u 


Arguing  as  we  did  for  the  functional  dependences  of  x and  u on  time, 
we  combine  Eqs.  ( 5 ) and  ( 2 ) to  obtain 


*k(t  + 


u f’'rk\  d 

V = -rk-ixk(t) + rk-i  Xk+1 


(t)- 


( 6 ) 


Equations  (4)  and  (6)  are  the  state  equations  for  an  arbitrary 
layer.  Because  of  boundary  conditions  at  the  surface  and  interface, 

we  must  derive  the  state  equations  for  X^(t+Tj)  and  xj^(t  + TK^  seParate^y* 
Seismic  source  signature,  m(t),  is  applied  at  the  surface.  In  this 
case  (compare  Figs.  6 and  4) 

U,  = x?  + ( 1 +r  ) m ( 7 ) 


Figure  6.  Components  of 

which  means,  from  Eq,  ( 1 ),  that 

Xj(t+Tl)  = rt  Xj(t)  + r^l  + rQ)  m(t)  ( 8 ) 

th 

For  the  K layer,  we  follow  customary  practice  to  assume  that 
no  signal  returns  up  from  the  basement;  hence,  (compare  Figs.  7 and  5) 


which  means,  from  Eq.  (2),  that 


Figure  7.  Components  of 


In  order  to  complete  the  state  space  model  description  for  the 
system  of  K layered  media  in  Fig.  1,  we  must  express  output  y(t)  in 
terms  of  input  m(t)  and  appropriate  states.  Output  y(t)  has  two  com- 
ponents (see  Figure  8),  a direct  reflection  of  m(t)  off  of  the  surface 
and  a transmission  of  an  upcoming  wave  through  the  surface.  That 
upcoming  wave  can  be  expressed  in  terms  of  a state  — x^  — as  depicted 
in  the  figure;  hence. 


Figure  8.  Components  of  y(t) 


( 12  ) 


y(t)  = rQ  m(t) 


B . Example 

In  order  to  illustrate  the  detailed  structure  of  state  equations  (11), 
we  write  them  for  the  two  layer  media  depicted  in  Figure  9 (K  = 2) : 

x^t  + Tj)  = r^^t)  + r j ( 1 + rQ)  m(t) 


d u 

Xj(t  + Tj)  = -rQ  x^t)  + r 


X2(t  + T2)  = 


-r0  xi(‘>  + ro  4 
r2  x2(t^  + r2  (-77)  4 


*2<*> 


4(t<-TZ>  = -T1  X2(t| 


The  output  equation  for  this  system  is  Eq.  (12).  For  the  purposes  of 
this  example,  we  have  chosen  to  order  the  states  by  layer.  Other 
orderings  are  possible  (see  Section  D).  Figure  9 depicts  primary  and 
multiple  reflections  and  illustrates  the  very  complicated  internal  be- 
havior of  the  system,  and  that  y(t)  has  components  which  occur  at 
specific  values  of  t.  For  compactness,  it  is  convenient  to  write  the 
four  equations  in  (13)  in  vector  matrix  notation,  as 


Xj(t+  Tj) 

x^(t  + Tj) 

X2(t  + T2) 
d, 

X2(t+T2) 


° r0(1  -r l>  7 r 1 


r2(l  + r1l/r1 


x“(t) 


ri(1  + ro) 


Figure  9.  Two  Layer  Example 


C.  A Compact  Description 

By  means  of  suitable  vector  matrix  notation,  it  is  possible  to  write 
Eqs.  (11)  and  (12)  in  a very  compact  manner.  To  begin,  we  define  the 
2K-dimensional  state  vector  x(t),  as 

x(t)  = col  [x^(t) , x^(t)  , x^(t) , x^(t) , . . . , x“(t) , x^(t)]  ( 15  ) 

Let  t_  be  the  K-dimensional  vector 

= col  ( Tj  , t2  , . . . , T ) ( 16  ) 

Then,  x(t+jT)  denotes  the  vector 

u d 

x(t  + jr)  = col  [x  (t  + t ) , x (t  + t ) , . . . , 


tk)  . x^(t  + Tv)]  ( 17  ) 


By  means  of  this  notation,  we  can  now  write  the  2K  state  equations 


in  (11),  as 


x(t  + T ) = Ax(t)  + B m(t) 


( 18  ) 


where  A and  B are  2K  X 2K  and  2K  X 1 matrices,  respectively.  Matrix 
A is  a band  matrix  with  zeros  along  its  main  diagonal;  its  structure 
for  K = 2 is  evident  from  Eq.  (14).  As  a further  example  of  its  structure 
we  give  A for  a three  layer  media  (K  = 3)  : 


r2(1  + rl)/ri 


0 -r 


0 r3(1+r2)/r2 


ro(1-ri)/ri 


0 ri(1-r2)/r2 


Matrix  B is  defined,  as 


( 19  ) 


B = col  f r x ( 1 + rQ)  , 0 , 0 , . . . , 0] 


Output  equation  (12)  can  now  be  written  in  terms  of  x(t)  , as 


where 


y(t)  = Cx(t)  + D m(t) 


C = ( 0 , -(l-r0)/r0  , 0 , , 0 ) 


( 21  ) 


( 22  ) 


D = r. 


( 23  ) 


It  is  important  to  understand  that  state  equation  (18)  is  a continuous - 
time  equation  that  involves  a vector  of  time  delays.  Some  remarks  on 
the  solution  of  such  an  equation  are  given  in  Section  D.  In  the  special 
case  of  equal  one  way  travel  times,  where  T =T  =...  = Tv  & t,  Eq.  (18) 

1 c.  IN. 


-li- 


Kii 


can  be  reduced  to  a vector  finite-difference  equation  by  choosing  t = kT: 

x[  (k  + 1)t  ) = A x(kT)  + B m(kT)  ( 24  ) 

In  the  more  general  case  of  nonequal  one  way  travel  times,  Eq.  (18) 
applies;  it  is  not  a finite -differ ence  equation,  but  belongs  to  the  class 
of  equations  known  as  causal  functional  equations. 

D.  Signal  Flow  Graphs 

In  order  to  develop  additional  insight  into  the  structure  of  our  state 
space  model,  it  is  instructive  to  portray  that  model  in  signal  flow  graph 

;Jc 

(SFG)  form.  Because  all  important  features  of  our  model  can  be 
observed  from  the  nature  of  a two  layer  example,  we  shall  develop  SFG's 
for  that  case,  leaving  generalizations  to  the  reader. 

The  state  equations  for  K=2  are  given  in  Eq.  (13),  and  the  obser- 
vation equation  is  given  in  Eq.  (12).  For  purposes  of  constructing  a 
SFG  it  is  useful  to  rewrite  Eq.  (13),  as  follows: 


*j(t) 


rixi(t ' V + ri(1  + V m(t ' Ti) 


( 25  ) 


Xl(t)  = -r0Xl(t-Tl)  + V(1-ri)/ri]X2(t-Tl) 
X2(t)  = r2X2(t_T2)  + r2[(1+  r i)/ril  xT(t  “ T2> 

X2(t)  = -rlX2(t-T2) 


A SFG  (which  is  an  interconnection  of  nodes  and  directed  branches,  with 

input  node  m(t),  output  node,  y(t),  and  intermediate  state  nodes  x'j'(t)  , 

d d 

Xj(t),  x^(t)  , x^t))  based  on  the  layered-ordering  of  the  states  in  Eq.  (25) 
is  depicted  in  Figure  10.  For  simplicity  in  drawing  the  SFG,  we  do 
not  show  the  explicit  dependence  of  node  signals  on  time;  and,  we  use 
operators  z}  and  z^  to  denote  Tj  and  sec.  time  delays,  respectively 


Signal  flow  graphs  are  essentially  the  same  as  block  diagrams,  but 
are  artistically  easier  to  draw.  See  Ref.  4 for  a discussion  on  SFG’s. 


-12- 


[i.e.,  z^m  = m(t-T  )].  Unfortunately,  this  SFG  does  not  have  a very 
revealing  structure,  and,  by  its  disorderly  appearance  suggests  that  a 
different  ordering  of  the  intermediate  state  nodes  might  be  more 
appropriate. 

A second,  and  much  more  revealing  SFG,  based  on  a computational 
ordering  of  the  states  in  Eq.  (25),  is  depicted  in  Figure  11.  We  refer 
to  the  order  of  the  states  - , x^  and  x^  - as  a computational 

ordering,  because  that  would  be  precisely  the  order  in  which  we  would 
solve  the  four  equations  in  (25).  We  will  have  more  to  say  about  the 
actual  solution  of  causal  functional  equations  in  a future  paper.  Observe 
that  the  Figure  11  SFG  has  a very  revealing  structure.  Three  feedback 
loops  are  present;  the  inner  loop  being  associated  with  layer  2,  the  other 
loops  being  associated  with  layers  1 and  2.  Another  example  of  a SFG 
based  on  a computational  ordering  of  states  is  given  in  Figure  12  for  a 
three  layer  media  (K=3).  This  SFG  illustrates  not  only  the  feedback 
paths  caused  by  layerings,  but  also  the  many  forward  paths  as  well. 

(i.e.,  there  are  6 feedback  paths  and  4 forward  paths),  and  serves  to 
point  up  the  very  complicated  nature  of  a layered  media  system. 

E.  Derivation  of  Reflection  Transfer  Function 

In  this  paragraph  we  develop  the  reflection  transfer  function  between 
input  rr(t)  and  output  y(t)  of  the  system  of  K layered  media  depicted  in 
Figure  1.  One  way  to  obtain  this  transfer  function  is  to  construct  the 
SFG  for  the  K layered  media  and  to  then  apply  Mason's  reduction  theorem 
(Ref.  4)  to  obtain  it  directly  from  the  topology  of  that  graph. 

A second  way  — the  one  which  we  shall  follow  — to  obtain  this 
transfer  function  is  to  develop  a system  of  recursive  equations  from  which 
it  can  be  calculated.  Our  approach  is  to  take  the  Laplace  transform  of 
the  2K  state  equations  (11)  and  the  output  equation  (12),  ultimately  obtaining 
the  Laplace  transfer  function  Y(s)/M(s).  We  shall  proceed  in  a systematic 
manner  focusing  our  attention  at  a layered-ordering  of  the  states. 
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a 


Let  X (s)  denote  a two-dimensional  Laplace  transformed  vector, 
k 


defined  for  the  kth  layer,  as 


X (s)  A col  [x£(s)  , xj*(s)  ] 


( 26  ) 


Additionally,  we  define  the  following  2x2  matrices 


Fk(s)  = 


eSTk 


C-SrO 


ST, 
e k 


( 27  ) 


( 28  ) 


Hk  = 


( 29  ) 


As  in  the  derivation  of  the  state  equations  (11),  we  must  separate 
the  analysis  of  the  first  and  last  layers  from  the  intermediate  layers 
when  we  take  the  Laplace  transform  of  Eq.  (11)  using  a layer-ordering 
of  states. 

For  k=l,  we  take  the  Laplace  transform  of  the  first  and  third 
(with  k=l)  equations  in  (11)  to  obtain 


Fj(s)  X^s)  = H j X2(s)  + 


( T"' ) 


( 30  ) 


For  k = 2,  3,  . . . , K- 1,  we  take  the  Laplace  transform  of  the  second  and 
third  equations  in  (11)  to  obtain 


Fk(s)Xk<s>  = Gk3k-1(a)  + HkXk+1<,) 


( 31  ) 


r 


Finally,  for  k = K we  take  the  Laplace  transform  of  the  second  (for  k = K) 
and  fourth  equations  in  (11)  to  obtain 


fk(s)  Vs)  = gk^k-i(s) 


( 32  ) 


Additionally,  from  Eqs,  (12)  and  (26),  we  obtain 


Y(s)  - rQ  M(s)  - 


(?) 


(0  1)  X,(s) 


( 33  ) 


from  which  we  see  that  to  compute  Y(s)/M(s)  we  must  first  determine 
Xj(s)  as  a function  of  M(s).  A recursive  procedure  for  accomplishing 
this  is  given  next. 

The  complete  set  of  Laplace  transformed  equations  for  the  system 
of  K layered  media  can  be  written,  as: 


fk  -k  gk  -k  - 1 


fk-i^k-i  = GK-l^K-2+HK-l^K 


FK-2  — K-2  = GK - 2 — K - 3 + HK - 2 — K - ! 


( 34a  ) 
( 34b  ) 
( 34c  ) 


f2x2  = g2x1  + h2x3 


Fl-Xl  = H1^2  + 


M 


( 34d  ) 


( 34e  ) 


where  the  argument  s has  been  omitted  for  notational  simplicity.  Observe, 

from  Eq.  (27),  that  ^(s)  is  invertible  for  all  k=l,2,  ...,K.  Solve  Eq. 

(34a)  for  X . substitute  the  resulting  expression  into  Eq.  (34b),  and 
K 

collect  common  terms  in  X , to  obtain  the  following  expression: 

K - 1 


(F 


K-l 


hk-ifkgk»2k-i  = gk-.2k-2 


( 35  ) 


Let  W.  denote  a 2x2  matrix,  defined  as 
J 


-17- 


t-  i 


W = F 
K K 


-1 


W.  = F.  - H.  W.  * G.  , 
J J J J+l  J+l 


( 36a 
( 36b 


where  j = K- 1 , K -2,  . . . , 1 . Then,  Eq.  (35)  can  be  rewritten  as 


WK- 1 *K-1  = GK-1  2K-2 


( 37 


Solve  Eq.  (37)  for  X , substitute  the  resulting  expression  into  Eq. 

(34c),  and  collect  common  terms  in  X _,  to  obtain  the  following 

K-  ^ 

expression: 


r-i 


(FK-2  ‘ HK-2  WK-1  ^-l*  -K-2  GK-2  ~K-3 


( 38  ) 


which  can  also  be  written,  as 


WK-2-K-2  GK-2-K-3 


( 39  ) 


Comparing  Eqs.  (39)  and  (37)  it  is  clear  that  we  can  continue  this  back- 
ward development  of  our  solution  for  Xj(s)  *n  same  manner.  For 
the  second  layer,  we  obtain 


W2  —2  = G2  — 1 


( 40  ) 


Solve  this  equation  for  X.^,  substitute  the  resulting  expression  into  Eq. 
(34e),  and  collect  common  terms  in  X^  , to  obtain 


(F 


-i  /ri<1  + ro)’\ 

1-  H1W2  G2>*1  = ^ q j 


M 


( 41  ) 


or. 


W 


( ri(1  + v \ 

( o ) 


M 


( 42  ) 


Finally,  solve  Eq.  (42)  for  Xj  and  substitute  the  resulting  expression 
into  Eq.  (33)  to  obtain  the  desired  reflection  transfer  function 


: 
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from  which  we  find 


2t?s 

e £ + r r 

+ 1 2 


1 2(t,+Tt)s  2tiS  2t9s 

e ' 1 2'  + r,r  e 1 + r r e 2 + r r 

12  10  0 2 


where 


2t,s  , 

..  ..  V 2 +r0r2/ri  ST, 

»i'2-  11  = iTT e 1 

e + r \r  2 


-aijl2.il  e“Tl 


( 47a  ) 


( 47b  ) 


_ _ 2tiS  _ 2t-)S 

Finally,  substitute  Eq.  (47)  into  Eq.  (43),  letting  z^=e  1 and  z^=e  ^ » 


to  obtain 


Y(s)  _ yZ_lZ2tr0^ir2ZltrlZZ  + r2 

“(■I  ' Vi'VzViVi'Vi 


( 48  ) 


which  is  the  desired  result.  For  the  special  case  where  Tj  = t2,  so 
that  z = z ~ z,  Eq.  (48)  simplifies  to 

1 w 


Y(s)  = V + (rltr0rlr2)Ztr2 

z2  + (r  r + r r )z  + r r 
'12  10  0 2 


( 49  ) 


which  is  precisely  the  same  result  derived  by  Robinson  in  Ref.  3. 


G.  Derivation  of  State  Transfer  Functions 

By  means  of  our  derivation  of  Y(s)/M(s),  it  is  also  possible  to 
obtain  X,(s)  as  a function  of  M(s).  Recall,  from  Paragraph  E,  that 

K 

wk*k=gk*k-i  ^ 


WK  - 1 — K - 1 = GK-1  — K-2 


W2  —2  G2  — 1 


W1  -1  = 


( 50  ) 
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Solve  the  last  equation  in  (50)  for  Xj*  substitute  the  result  into  the 
next  to  the  last  equation,  solve  it  for  * efc*  t°  obtain  the  following 
state  transfer  function 


1 i 

Xk(s)  =iSk  [wi  (8)  M(8) 


where 


So,  for  example,  when  K = 2, 


and 


X = Wj  Gj  M 


X2  = Wz‘  G2  w;‘  G1  M 


( 51  ) 


( 52  ) 


( 53a  ) 
( 53b  ) 


By  means  of  the  state  transfer  functions,  it  is  possible  to  compute 
the  states  of  a K layered  media  directly  in  terms  of  the  seismic  input 
signature  m(t). 


III.  RELATIONSHIP  OF  NEW  COMPLETE  MODEL  TO  TRADITIONAL 

MODEL 

A.  Introduction 

In  this  section  we  develop  a state  equation  representation  for  the 
traditionally-derived  model  of  a system  of  K layered  media  (Ref.  3,  for 
example);  however,  unlike  the  more  traditional  derivations,  we  shall 
assume  unequal  one  way  travel  times.  Then,  by  means  of  a two-layer 
example,  we  establish  the  coordinate  transformation  matrix  which  relates 
the  traditional  and  new  models. 

B.  Derivation  of  a Traditionally-Derived  State  Space  Model 


The  starting  point  for  the  development  of  the  traditional  model  is 
the  Figure  13  ray  diagram  (taken  from  Ref.  3).  We  use  the  (traditional) 


Figure  13.  Reflected  and  transmitted  waves  at  interface  k 


4 


symbols  and  d^  to  denote  the  upgoing  and  downgoing  waves  in  the 
kth  layer,  respectively;  and,  adopt  the  convention  that  waves  at  the  top 
of  a layer  occur  at  present  time,  t.  To  develop  the  traditional  model 
we  direct  our  attention  at  the  intersection  point  of  the  ray  diagram  and 
apply  superposition  to  obtain  the  following  equations  for  signals  u^(t  + t^) 
and  d^^ft),  which  leave  that  point; 


\lt  + V = rkdk(t-Tkl  + <1-rk)uk+i(t)  1 1 

dk+i(t)  3 ,ltrk' dk(t-V  - pk"k+i|t)  ( 2 

These  equations  are  applicable  for  k = 1 , 2, . . . , K- 1 . At  the  surface 
(Figure  14a),  we  obtain 


u0(t)  = rodo(t)  + (1'ro,ui(t)  ( 3 ) 

dj(t)  = (1  + rQ)  dQ(t)  - rQ  u^t)  ; ( 4 ) 
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( a ) . ( b ) 

Figure  14.  Reflected  and  transmitted  waves  at  ( a ) surface  and  (b ) interface  K. 


and,  at  the  K**1  interface,  we  assume  that  u^+j(t)  = 0,  to  obtain  (Figure 
14b) 


lK|t  + V = 'kV'V 

( 5 ) 

W’  = (1+rK,dK(t-’K» 

( 6 ) 

Signal  Ug(t)  in  Eq.  (3)  is  the  measurable  system  output;  i.  e.  , u^ft) 

is  analogous  to  y(t)  in  Section  II.  Signal  d (t)  is  also  a system  output; 

K+ 1 

but,  since  it  cannot  be  measured,  we  shall  ignore  it  in  following  analyses. 


Signal  d^(t),  the  seismic  source  signature,  is  analogous  to  m(t)  in  Section 


II. 


It  is  convenient  to  group  Eqs.  (1),  (2),  (4)  and  (5)  in  a layer  ordering, 


as  follows  : 


djft)  = -roU^  + O + V V*> 


u^t  + Tj)  = Tjd^t  - Tj)  + (1  - rL)  u2(t) 


d.(.)  = (i  + q.j)  -Tj.,)  - rj.i  uj<" 


u.(t+T.)  = r.  d.(t-T.)  + (1  -r.)  u.  (t) 
J J i J J J J+l 


j = 2.3, 


K-l 


The  reader  may  wonder  why  we  list  the  downgoing  signal  in  each  layer 
first,  whereas  in  Section  II  we  listed  the  upgoing  signal  first.  The  reason 
for  doing  this  will  become  apparant  in  Paragraph  C. 


K 


K-l  K-l  K-l  K-l  K' 


Vt+,K>  = 'kV'-’k1 


( ^ ) 


This  system  of  2K  equations  is  not  in  state  equation  format,  yet,  since 

signals  in  its  left-hand  side  occur  at  t and  delayed  times,  and  signals 

on  the  right-hand  side  occur  at  t,  t - t.  , and  t - t..  In  order  to  put 

J-l  J 

Eq.  (7)  into  state  equation  format,  let 


e.(t  + t.)  = d.(t) 
J J 1 


( 8 ) 


for  all  j = 1,  2,  . . . , K.  Then,  Eq.  (7)  becomes 


ei(t  + Ti)  = Vi(t,  + (1  + ro)do(t) 


U1  (t  + Ti>  = + (1  -r})  u2(t) 


e.(t+T.)  = d + r.^)  e.Jt)  - r.^  u.(t) 

u.(t  + T.)  = r.e.(t)  + (1  - r.)  u (t) 

J J J J J J+l 

eK(t  + V = ,1+rK-l’  'K-l**'  - rK-l  UKll) 


j = 2,  3, 


, K- 1 


UK(t+V  = rKeK(t) 


( 9 ) 


By  means  of  transformation  (8)  each  pair  of  equations  in  (9)  now  only 

involves  two  time  points,  t+T.  and  t. 

J 

Let 

X(t)  = col  [e^t)  , u^t)  , e2(t)  , u2(t) eR(t)  , uR(t)]  ( 10  ) 

and 

^(t  + T)  = col  [e  (t  + Tt)  , u^t  + Tj) eK(t  + TK>  ' UK(t+  tK^  ( 11  ) 


Then,  Eq.  (9)  can  be  written  as 

*(t  + T)  3 Ax  X(t)  + dQ(t) 


( 12  ) 


where  A and  B are  2K  X 2K  and  2K  X 1 matrices,  respectively.  Matrix 


A is  a band  matrix  with  zeros  along  its  main  diagonal.  It  has  zero  and 

X 

non-zero  entries  in  exactly  the  same  locations  as  does  matrix  A [see  Eqs. 


(11-14)  and  (11-19)  for  examples].  Matrix  B is  defined  as 

x 

= col  [1  + r^  , 0,  0 0] 


( 13 


Additionally,  output  equation  (3)  can  be  written  in  terms  of  X_(t),  as 

( 14  ) 


where 


and 


V*’  ‘ cx**,  + Vo*" 


(0,1-r  ,0 0) 

0 


D = r 
X o 


( 15  ) 


( 16 


In  crder  to  illustrate  the  detailed  structure  of  Eq.  (12),  and  especially 
the  structure  of  A , we  write  out  that  equation  for  a two  layer  media: 


0 

-r 

0 

0 

° \ 

/ ei(t) 

rl 

0 

0 

'"r‘  1 

[ Ul(t) 

1+rl 

0 

0 

r' 

1 e2(t) 

0 

0 

r2 

0 J 

\ U2(t) 

\ 


e I (t  + T j ) 
Uj(t  + T j ) 

e2(t+T2) 

U2(t+T2} 


( 17  ) 

C.  Relationship  of  Models 

We  have  derived  two  state  space  models  for  the  same  linear  system: 


and 


N 


x(t  + t)  = A x(t)  + B m(t) 
y(t)  = C x(t)  + D m(t) 


X.(t  + D = Ax>L(t)  + B^  dQ(t) 

V*>  = cx  1(t)  + Dx  do(t) 


( 11-18  ) 

( H-21  ) 

( III- 12  ) 
( III- 14  ) 


in  which  m(t)  = d^(t)  and  y(t)  = u^(t).  Because  of  the  linear  nature  of 
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and  S^,,  x(t)  and  ^_(t)  must  be  related  by  a linear  transformation,  say  T, 


c « » 


x(t)  - T£(t) 


( 18  ) 


We  wish  to  determine  the  exact  structure  of  T so  that  we  can  make  a 
physical  connection  between  state  vectors  x(t)  and  ^(t). 

If  Eq.  (18)  is  true,  then  and  S,^  are  equivalent  if 


A = T A T" 

X 

B = T B 

X 


C = C T 
X 


( 19  ) 
( 20  ) 


( 21  ) 


We  have  before  us  the  following  problem:  given  matrices  A,  A , B,  B , 

X X 

C and  C , find  matrix  T such  that  Eqs.  (19),  (20)  and  (21)  are  simul- 

X 

taneously  satisfied.  A direct  approach  to  solution  of  this  problem  is  to 

rewrite  (19)  and  (20)  as  AT  = TA^  and  CT  = C^,  and  solve  these  equations 

plus  Eq.  (20)  for  the  elements  of  T.  This  is  a large  system  of  linear 

equations;  but,  it  can  be  solved  in  this  direct  manner.  A second,  and 

much  easier  approach  to  solution  of  this  problem  is  to  compare  the 

structures  of  matrices  A and  A , B and  B , and  C and  C , to  observe 

XX  X 

that  nonzero  elements  in  comparable  matrices  occur  in  exactly  the  same 
location;  hence,  T must  be  a diagonal  matrix. 


Consider  a two  layer  media  system,  where  A,  B and  C are  given 

in  Eqs.  (11-14)  and  (11-12),  and,  A , B and  C are  given  in  Eqs.  (17) 

XX  X 

and  (15).  Let 


T = diag  (tn  , t22,  t33,  t44) 


( 22  ) 


then,  from  the  equation  B = TB^,  we  find  the  equation  CT  = C , 

we  find  t _ = - r • and,  from  AT  = TA  we  find  t_ ,=  r,  and  t = -r  ; hence, 

22  0 X 33  2 44  1 


= dia8  <V  *r0  ' r2'  -V 


( 23  ) 


Using  the  definitions  of  x(t)  and  ^.(t),  and,  Eq.  (23),  we  establish, 
from  Eq.  (18),  that 
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xj(t)  = = r1d1(t  ’ V 


xl(t)  = -rgU^t) 


X2(t)  = r2C2(t)  = r2d2(t‘T2) 


x2(t)  = -riU2(t) 


( 24  ) 


By  means  of  these  simple  relationships,  we  relate  the  new  and  traditional 
states  as  dipicted  in  Figure  15.  Observe  that  the  states  in  the  model 
derived  in  Section  II  are  the  reflections  of  the  states  in  the  model  of 
Section  III.  In  fact,  if  we  compare  the  following  versions  of  Eqs.  (II- 1) 
and  (II-2)  with  Eq.  (24), 


Xl(t)  = rlUl(t  ‘ Tl) 


xi(t)  = •ro"I(t-Ti) 


X2(t)  = r2U2(t'T2) 


x2(t)  = -t2) 


( 25  ) 


we  see  that  d =u  , d -yi  , u = 5.  and  u = % , which  further  supports  this 

11  Cbli  C*  L. 

conclusion. 

Generalizations  of  the  preceding  results  are  obtained  in  a straight- 
forward manner.  For  a K layer  media, 

T = diag  (r  1 , -rQ  , r^,  rR , -rK1)  • ( 26  ) 


IV.  A PRIMARY  REFLECTION  MODEL 


A.  Introduction 

In  certain  applications  it  is  useful  to  approximate  the  complete  seis- 
mogram by  just  the  primary  reflections.  Conceptually,  the  primary  re- 
flections (i.e.,  primaries)  are  easily  obtained  directly  from  ray  diagrams; 
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Figure  15.  Relationships  between  signals  in  the  two 
models. 

however,  that  approach  is  not  very  useful  in  estimation  procedures 
where  one  would  be  interested  in  associating  a dynamical  model  with 
the  primaries.  In  this  section  we  demonstrate  how  to  modify  the 
complete  model  in  a very  simple  manner  so  as  to  obtain  a state  space 
model  which  generates  only  the  primaries.  We  refer  to  this  model 
either  as  a primary  reflection  model  or  as  a primaries  model. 

B.  Derivation  of  Primaries  Model 

State  equations  for  a primaries  model  are  obtained  from  Eq.  (II- 11) 

by  deleting  the  term  r x^(t)  in  the  equation  for  x^(t  + T ),  and  the  term 

d ^ ^ ci  A 

r x (t)  in  the  equations  for  x (t  + t ),  k = 2,  3,  . . . , K.  Our  reason  for 
k k K K 

deleting  these  terms  is  based  on  the  following  observation:  a multiple 

reflection  can  occur  only  when  an  upgoing  wave  reflects  off  of  the  top 

d 

of  a layer.  Observe,  in  Figure  4 that  the  component  x^  of  represents 

just  such  a reflection,  as  does  the  component  x^  of  ^ (Figure  6);  hence, 
these  components  are  deleted  from  the  complete  model  in  our  primaries 
model,  whose  equations  are  given  below. 
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r 


X1 1 (t  + Ti>  = r j ( 1 + r0)  m(t) 


rk[<1+rk-il/rk-i]xi,k-i1*1 


k = 2,  3,  . . . , K 

^k(t+Tk» = -tk-ixikl*)  + rk-i[|1-rk)/rklx?.k+i(t) 


( 1 ) 


k = 1,  2,  . . . ,K-1 
X?K(t+TK)  = -rK-lXlK(t) 


yj(t)  = r0  m(t)  -[(i-rQj/'olXjjft) 


( 2 ) 


u d 

In  these  equations  x , x and  denote  primary  upgoing  and  downgoing 
states  and  the  primary  reflection  portion  of  the  complete  output,  respectively. 


C.  SFG  for  Primaries  Model 

An  important  property  of  the  primaries  model  is  that  its  SFG  has 
no  feedback  paths.  We  illustrate  this  for  a three  layer  media  in  Figure 
16.  This  property  causes  y j (t ) to  be  comprised  of  a finite  number  of 
terms  — K+l,  for  a K layered  media,  and  enables  y ^ (t ) to  be  computed 
in  a very  straightforward  manner. 

D.  Example 

In  order  to  validate  our  primaries  model,  we  consider  a three  layer 
media  and  compute  y ^ (t ) in  two  different  ways:  ray  tracing  and  Eqs.  (1) 

and  (2).  A ray  tracing  solution  is  depicted  in  Figure  17.  We  assume 
that  the  reader  is  familiar  with  the  details  of  such  a solution;  hence,  we 
do  not  label  the  amplitudes  of  the  rays  on  that  figure.  Signal  y ^ (t)  has 
four  components  a,  b,  c and  d which  appear  in  that  order  in  the  following 
expression  for  y ^ (t ): 

y (t)  = rn  m(t)  + r (1  - r^)  m(t  - 2t  ) 


+ r2(1  ' rQ*  (1  “ rl)  m<t-2T1-ZT2> 


■l: 


Figure  17.  Ray  diagram  depicting  primary  reflections  as 
transmitted  to  the  surface. 


+ r3(l  - r^)  (1  - rj)  (1  - r^)  m(t  - 2-^-  2t  Z-r^)  ( 3 ) 

The  solution  for  y^(t)  from  Eqs.  (1)  and  (2)  can  be  obtained  either 
directly  from  the  Figure  16  SFG  or  from  those  equations.  Because  we 
shall  need  some  of  the  states  for  future  calculations  in  Section  V,  we 
shall  list  the  states  obtained  in  order  of  solution  from  Eq.  (1): 


.-,4 

x^(t  + Ti)  = r j(l  + rQ  ) m(t)  ( 4a  ) 

X12(t+  T2>  = r2(1  + r0)  (1  + rl)  m(t  " Tl)  ( 4t>  * 

x^3(t  + t3)  = r3(l  + r2)  (1  + rL)  (1  + rQ)  m(t  -Tj-T2)  ( 4c  ) 

Xn(t  + T3)  = 'r2r3(1  +r2)(l  + ri)(l  +r0)  m(t-  Tj-t2-t3)  ( 4d  ) 

X?2(t  + T2)  = -rir2(1+rl)(1  + r0)  m(t-TrT2> 

2 

‘rir3<1  * r2)^1  + ri^1  + ro)  "Tl'T2"  2t3)  ( 4e  ) 
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Xll(t+V  = -r0ri(1+r0)  - r0r2(1  'ri)(1+r0)  m(t-Tr2T2) 


- rQr3(l  - r2)  (1  - r2)  (1  + rQ)  m(t  - Tj- 2t2- 2t3)  ( 4f  ) 


•Substitute  x^ft)  into  Eq.  (2)  to  obtain  precisely  the  same  expression  for 
y ^ (t)  as  obtained  in  Eq.  (3)  via  ray  tracing. 

E.  Primary  Reflection  Transfer  Function 

In  Section  II-E  we  derived  a recursive  procedure  for  computing  the 
transfer  function  Y(s)/M(s).  In  general,  y(t)  will  contain  an  infinite 
number  of  terms  — primaries  and  all  multiples.  This  is  due  to  the  fact 
that  the  complete  transfer  function,  Y(s)/M(s),  is  a ratio  of  two  poly- 
nomials [see  Eqs.  (11-48)  and  (11-49),  for  examples].  The  primaries 
transfer  function,  Y^(s)/M(s),  on  the  other  hand  has  no  poles.  This 
is  evident,  for  example,  from  the  example  in  Figure  16. 

We  do  not  wish  to  redevelop  the  recursive  procedure  for  computing 
Y j (s)/M(s)  in  the  detail  of  Section  II-E.  The  results  will  look  exactly 
like  those  given  in  Eqs.  (11-43)  and  (11-36),  except  that  all  matrices  and 

vectors  should  have  an  additional  subscript  1 (e.  g.  , W , Y ) to  associate 

1 lx  1 

them  with  the  primaries  model.  An  important  difference  occurs  in  the 

definition  of  F,,(s): 
lk 


Matrices  F„  (s)  are  all  lower  triangular;  and, 
lk 

strate,  by  direct  calculations,  that  matrices  W 
W are  also  lower  triangular;  hence,  is 


one  can  therefore 
W 

l.K-1’  " 1 , K-2  ’ *' 
lower  triangular. 


demon - 
. . W. 


12 


say 


W 


-1 
1 1 


(s) 


a(s) 

0(8) 


( 6 ) 


F’rom  Eqs.  (11-43)  and  (6),  we  see  that 
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which  proves  our  assertion  that  the  primaries  transfer  function  has  no 
poles.  We  shall  have  more  to  say  about  this  property  of  the  primaries 
model  in  Section  VI. 

V.  A SECONDARY  REFLECTION  MODEL 


A.  Introduction 

In  Section  IV  we  developed  a primaries  model  which  permits  us  to 
remove  the  primary  reflections  from  the  complete  output.  We  are 
interested  in  learning  whether  it  is  possible  to  develop  a secondary 
reflection  model  (or,  secondaries  model)  which  would  generate  just  the 
secondary  reflections.  Secondary  reflections,  which  are  the  first  multiple 
reflections,  are  components  of  the  complete  output  which  are  due  to 
exactly  three  reflections  within  a K layer  media  system.  In  some 
applications  the  secondaries  can  be  significant  contributors  to  the  complete 
output;  whereas,  in  other  applications  they  may  be  of  such  small  amplitude 
(due  to  the  product  of  three  reflection  coefficients  in  each  secondary)  that 
they  can  be  neglected.  In  both  situations  it  is  of  interest  to  establish 
whether  they  should  or  should  not  be  accounted  for.  The  secondaries 
model  which  we  develop  in  this  section  will  permit  one  to  focus  his 
attention  just  at  the  secondary  reflections;  hence,  he  will  be  able  to 
establish  the  significance  of  the  secondaries  without  the  encumbrance  of 
primaries  and  other  multiples,  which  are  present  in  the  complete  output. 

Let  x^.(t)  and  x^.(t)  denote  secondary  upgoing  and  downgoing  states 
th  ^ ** 

in  the  j layer  (j  = 1 , 2, . . . , K),  and,  let  y ^(t)  denote  the  secondary  re- 
flection portion  of  the  complete  output.  Our  objective  in  this  section  is 

u d 

to  develop  state  equations  for  x^j(t)  and  x^(t)  and  the  equation  for  y^t). 
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B.  A First  Residual  Model 


To  begin,  we  form  a residual  model  between  our  complete  model 
and  primaries  model.  Subtract  the  respective  equations  in  Eqs.  (II- 1 1 ) 
and  (IV- 1)  to  obtain: 

Xj(t  + Tj)  ~x^1(t  + T1)  = r!  xj(t)  ( la  ) 

Xk|t  + V-Xik(t+Tk)  = rkxk(t)  + rkf(1+rk-i)/rk-i1 


■[<_,<*>  - 


k-r  ' i,k-r 


( lb  ) 


k = 2,  3,  . . . , K 


*k(t  + V - xik,,tTk>  = -rk-i[*k,t|  • bk"" 


+ rk-i[<1-rk,/rk![xk+i<‘»-  4.k+i<,,!  ( lc  1 


k = 1,2,..  . ,K-1 


4(‘  + V • = -rK-ltXK<t»  • XlK<t)I 


Additionally,  subtract  y^(t)  [Eq.  (IV-2)]  from  y(t)  [Eq.  (11-12)],  to  obtain 
y(t)  - y J (t)  = -[(1  -r0)/rQ]  [ x^(t)  - x^(t)  ] ( 2 ) 

u d 

Next,  we  define  first  residual  states  and  anc^  first 

residual  output,  Sj(t),  as 


= x^t)  - x«(t) 
§jj(t)  = xf(t)  - x^.(t) 


Sj(t)  = y(t)  - y j (t) 


( 3a  ) 


( 3b  ) 


d d 

Adding  and  subtracting  the  term  r x (t)  from  Eq.  (la)  and  r x (t)  from 

111  K IK 

Eq.  (lb),  and  substituting  Eqs.  (3a)  and  (3b)  into  Eqs.  (la)-(ld),  we  obtain 
the  following  first  residual  state  equations: 


5U<,  + Ti)  = ri?n(t)  + ri*n(t) 


!ik(,  + Tk>  = rk5ik(t|  + rkt(1+rk-i,/rk-il!l‘,k-i<t)  + Vik"1 

k = 2,  3,  . . . , K 

!u(t+V  = •rk-i!ik<t)  * rk-it(1-rk)/rklpn,k«(t) 

?lK<t  + TK’  = -rK-l  ?“K(t) 

Additionally,  the  first  residual  output,  s^(t),  is 

«!<*)  ?U(t) 


( 5 ) 


( 6 ) 


At  this  point  some  important  observations  are  in  order.  By  its 
very  definition,  s^(t)  comprises  all  the  multiple  reflections  for  a K 
layer  media  system.  The  first  residual  model  (5)  and  (6)  can  there- 
fore be  used  to  generate  the  complete  set  of  multiple  reflections  — 
secondaries,  tertiaries,  etc.  Observe  also  that  the  first  residual  state 
equations  are  driven  by  the  dawngoing  states  from  the  primaries  model; 
hence,  the  primaries  model  is  coupled  into  the  first  residuals  model, 
but  in  an  open-loop  manner. 


C.  Secondaries  Model 


Compare  state  equations  (5)  with  state  equations  (II- 11),  to  observe 


that  they  are  structurally  quite  similar.  Recall  that  in  the  development 


of  the  primaries  model  we  deleted  the  terms  r x^(t)  and  r x^(t)  from  the 

11  K K 

U 

state  equations.  We  conjecture,  therefore,  that  the  state  equations  for 


the  secondaries  model  are  obtained  from  Eq.  (5)  by  deleting  the  term 
r 1 ?<1 1 (t)  *n  equation  for  £^j(t+T^),  and  the  term  *n  the  equations 

,K.  The  resulting  secondaries  model  is: 


f°r  ?lk(t  + Tk)'  k=2’3’* 


x^(t  + T,)  = r,x?  (t) 


21' 

u 


1 11 


Vt+V  ■ rJ|1+r 


k-l)/rk-l1*2.k-l,,)  + rkXlk(,) 


1 


. 

i 


*tt(t  + V = -rk-lX2k(t)  + rk-l[‘1-rk)/rb,Xlk+l(tl 


( 7 ) 


k = 1,2,..  . ,K-1 

X2K(t+V  = -rK-lX2K(t) 


and 


y2(t)  =.((l-r0)/r0l  x21(t) 


( 8 | 


A proof  of  this  result  will  appear  in  a forthcoming  paper.  We  validate 
it  by  means  of  an  example  in  Paragraph  E.  We  summarize  the  essence 
of  the  secondaries  model  in  Figure  18.  Observe,  again,  the  dependence 
of  the  secondaries  model  on  all  of  the  downgoing  states  from  the  pri- 
maries model.  Additional  discussions  on  the  structure  of  the  block 
diagram  in  Figure  18  are  given  in  Section  VI. 


D.  SFG  for  Secondaries  Model 

An  important  property  of  the  secondaries  model  is  that  its  SFG  has 
no  feedback  paths.  We  illustrate  this  for  a three  layer  media  in  Figure 
19.  This  property  causes  y ^(t)  to  be  comprised  of  a finite  number  of 
terms  and  enables  y ^(t)  to  be  computed  in  a very  straightforward  manner. 


E.  Example 

To  validate  our  secondaries  model,  we  consider  a three  layer  media 
and  compute  y ^(t)  in  two  different  ways:  ray  tracing  and  Eqs.  (7)  and  (8). 


£ 

m(t) 

Primaries 

Secondaries 

] 

Model 

Model 

I 


y,(t) 


y2(t) 


Figure  18. 


Generation  of  primaries,  y (t), 

y2(t)- 


and  secondaries. 
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The  ray  tracing  solution  is  much  more  complicated  for  the  secondaries 
than  for  the  primaries  (Section  IV-C);  but,  can  be  accomplished  for  the 
three  layer  media  by  superposition  as  demonstrated  by  the  six  situations 
depicted  in  Figure  20.  Signal  y2(t)  has  14  components  a,  b,  c,  . . . , and  n 

which  appear  in  that  order  in  the  following  expression  for  y?(t):  let 

2 2 2 
a = 1 - rQ  , 0 = 1 - r and  y = 1 - r^  , then 

2 

y 2(t)  = - rQ  r j a m(t-4T])  - rQr } r 3 m(t  - 4^-  2t2> 

- ror1r3a0Y  m(t  “4Tr  2t2'  2t3) 

2 2 

- rQTjT^a  p m(t- 41^- 2t2)  - r^r^S  m(t-4T1-4T2> 

2 

- r0r3r2aP  Y m(t  - 4tj  - 4t2-2t3> 

- r^r  j r^a  0 y m(t  -4t1_2t2-2t3) 

2 

- r0r2r3a^  y m(t  - 4t^- 4t2-  2t3> 

2 2 2 

- r^ r^a  3 y m(t  - 4^- 4t2- 4t3) 

2 

- r^a  p m(t  - 2t1*4t2> 

- r}  r2r3a  3 y m(t  - 2tj  - 4t2- 2t3) 

- r j r 2r  3a3  y m(t  - 2t1-4t2-  2t3) 

- Tjr3a3  y2  m(t  - 2Tj  - 4t2- 4t3> 

- r2r3«  3 V m(t  *2ti‘2t2'4t3)  (9) 

Many  of  these  terms  can  be  grouped  together;  but,  as  displayed  we  are 
able  to  correlate  each  term  with  the  secondary  reflection  components 
depicted  in  the  six  figures  of  Figure  20. 

The  solution  for  y 2 ( t ) from  Eqs.  (7)  and  (8)  can  be  obtained  directly 
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Figure  20.  Ray  diagrams  depicting  secondary  reflections  as  trans- 
mitted to  the  surface.  Dashed  rays  indicate  portion  of 
a primary  reflection  which  generates  subsequent 
secondary  reflections. 


from  those  equations,  or  from  the  SFG  in  Figure  19.  We  give  some 
of  the  details  for  this  latter  approach  next.  Treating  zj  » z2  anc*  zi  as 


delay  operators,  we  establish  that 

y2(t)  ~ -C(1-ro>/ro]  X21^ 


( 10  ) 


X21(t)  = t-r0Vrr0r2MlVr0r3evZlZ2Z3]xn(t) 

+ t-r0r2(1-ri)ZlZ2  -r0r3(1-ri>VzlZ2Z3^XU(t) 
+ [-r0r3(l-r1)(l-r2)ziz2z^]x^3(t) 


( 11  ) 


where,  from  Eqs.  (IV -4f ) , (IV-4e)  and  (IV-4d), 

d 2 2 2 

Xll(t)  = t-r0ri(1+r0)zi  - r0r2(1+r0)&zlz2 

2 2 2 

'r0r3<1  +r0)PvZlZ2Z3  1 m(t) 


( 12a  ) 


X12(t)  = f‘rir2(1  + r0)(1  + ri,ZlZ2‘ rir3(1  + r0,(1  + rl)yZlZ2Z3^m(t) 


( 12b  ) 


Xn(t)  = [-r2r 3(  1 + rQ )( 1 + TjKl  + r2)  ZjZ2z3  ] m(t) 


( 12c  ) 


Substitute  Eqs.  (12)  into  Eq.  (II)  and  carry  out  the  resulting  multiplications 
to  obtain  y2(t)  exactly  as  given  in  Eq.  (9).  The  first  term  in  Eq.  (11) 
yields  the  first  9 terms  of  y2(t);  t^ie  sec°nd  term  in  Eq.  (11)  yields  the 
next  4 terms  of  y t ) ; and,  the  last  term  in  Eq.  (11)  yields  the  14th 
term  of  y 2<t ) . 


F.  Remark 

It  is  only  in  simple  examples  which  involve  a small  number  of  layers 
that  one  can  obtain  y ^(t ) via  ray  tracing  techniques.  There  are  just  too 
many  possibilities  for  enumeration  in  larger  systems.  However,  our 
secondaries  state  space  model  is  applicable  regardless  of  how  many 
layers  are  present.  We  are  presently  studying  computer  implementations 
of  these  equations. 


VI.  A COMPLETE  DECOMPOSITION  OF  A SEISMOGRAM 


) . 


f 


L.  . 


Our  results  in  the  two  preceding  sections  suggest  that  a complete 
decomposition  of  a seismogram  [i.  e.  , y(t)],  into  a superposition  of 
outputs  from  primaries,  secondaries,  tertiaries,  etc.  models  is  possible. 
We  will  not  go  into  all  of  the  details  of  such  a decomposition  here;  but, 
instead  shall  sketch  its  development. 

To  begin,  let  us  sketch  the  development  of  a tertiaries  model. 
Tertiary  reflections  are  the  second  collection  of  multiple  reflections, 
and  are  those  components  of  the  complete  output  which  are  due  to  ex- 
actly five  reflections  within  a K layer  media  system.  To  obtain  a 
tertiaries  model,  we  proceed  exactly  as  in  our  development  of  the 
secondaries  model.  First,  we  form  a second  residual  model  between 
the  first  residual  model  in  Eqs.  (V-5)  and  ( V - 6 ) and  the  secondaries 


model  in  Eqs.  (V-7)  and  (V-8).  The  second  residual  state  variables 


are  ?^.(t)  and  ?_.(t),  where 
2J  2J 


?2j(t>  = ?lj(t)  - *2^ 


( la  ) 


and 


4^  - ?ij(t)  - x2j{t) 


( lb  ) 


for  j = l,2, . . . , K,  and,  the  second  residual  output  is  s2(t),  where 


s2(f)  = Sjft)  - y2(t) 


( 2 ) 


The  second  residual  state  equations  will  look  exactly  like  the  first 


residual  state  equations  in  (V-5)  with  an^  replaced  by  % and 

i j d ^ ^ ^ ^ 

5^.,  respectively,  and  x^  ^ and  x^  replaced  by  x^^  and  x^*  respectively. 


The  second  residual  state  equations  are  driven  by  the  downgoing 
components  of  the  secondaries  state  vector. 

The  tertiaries  model  is  obtained  from  the  second  residual  model 


by  deleting  the  terms  r^x^  (t)  and  r^x^(t)  in  the  ?2j(t  + Tp  an<*  + Tk^ 


state  equations.  The  resulting  tertiaries  model  is: 
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J 


= ri  X21(t) 


X3k(t+Tk)  rkf(1  + rk-l)/rk-l^X3,k-l(t)  + rk  X2k 


(t) 


and 


k = 2,  3,  ....  K 

X3k(,tTk'  = -rk-l  X3k^  + rk-l^(Irk^/rk^X3,  k+1^*^ 
X3K,,+  ,K)  = -rK-lX3K,') 


y3«)  = -t(l-r0)/r0Ix3](t) 


V ( 3 I 


( 4 ) 


Observe  that  the  tertiaries  model  is  driven  by  the  downgoing  components 

of  the  secondaries  state  vector.  Figure  18  can  now  be  expanded  to 

d 

include  the  tertiaries  model  driven  by  x^(t)  with  output  y^(t). 

In  this  manner  we  can  form  successive  residual  models  and 

subsequent  n~aries  models,  each  model  being  of  dimension  2K.  The 

following  complete  (canonical)  decomposition  of  y(t)  has  been  demonstrated 

by  the  second  author  of  this  paper: 

The  complete  output,  y(t),  from  a K layer  media  system, 
which  is  comprised  of  the  superposition  of  primaries, 
secondaries,  tertiaries,  etc.  , can  be  obtained  from  a 
single  model  of  order  2K  — the  complete  model  — or 
from  an  infinite  number  of  models,  each  of  order  2K, 
interconnected  as  shown  in  Figure  21. 

There  are  many  interesting  implications  and  potentially  useful 

applications  for  this  decomposition.  In  seismic  applications,  where 

reflection  coefficients  are  often  quite  small  [|r^|<0.3],  the  decomposition 

could  be  used  as  a computational  tool  for  approximating  the  complete 

seismogram  y(t)  by  a small  collection  of  constituent  reflections,  such  as 

y (t)  + y (t).  Each  one  of  the  constituent  models  in  Figure  21  is  a non- 
1 ^ 

feedback  model,  whereas  the  complete  model  is  a very  complicated 
feedback  system;  hence,  computation  via  a truncated  decomposition  of 
y(t)  may  be  quite  expeditious.  Additionally,  by  means  of  the  canonical 
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decomposition  it  should  be  easy  to  compare  the  relative  importance  of 
secondaries  to  primaries,  etc.  It  should  also  be  possible  to  use  the 
primaries  model  (or  a combination  of  primaries  and  secondaries  models) 
in  deconvolution  techniques,  which  often  assume  that  primaries  are  the 
only  significant  terms  in  y(t).  Finally,  it  should  be  possible  to  depict 
primaries,  secondaries,  etc.  in  different  colors  on  a synthetic  seismo- 
gram, when  the  seismogram  is  constructed  via  the  canonical  decomposition. 

VII.  NON-NORMAL  INCIDENCE  MODELS 

All  of  the  preceding  developments  were  for  the  situation  of  normal 
incidence.  In  this  section  we  suggest  a way  in  which  these  results  can 
be  extended  to  the  situation  of  non-normal  incidence.  Our  interest  in 
obtaining  such  an  extension  is  motivated  by  the  arrays  of  geophones  or 
hydrophones  which  simultaneously  record  seismograms  in  response  to  a 
seismic  source  signature.  Except  for  a sensor  located  exactly  at  the 
source  point,  the  other  sensors  receive  non-normally  incident  signals. 

We  wish  to  extend  our  state  space  models  in  as  simple  a manner  as 
possible  to  describe  these  other  sensors.  Our  approach  is  to  obtain 
an  expression  for  a travel  time,  (i  = 1 , 2, . . . , K),  which  we  refer  to 
as  one-way  non-normal  incidence  travel  time  for  horizontal  offset 
distance  d.  This  travel  time  can  then  be  use  in  place  of  the  one-way 
normal  incidence  travel  time,  t.  (i  = 1 , 2, . . . , K),  in  all  of  our  preceding 
models. 

For  normal  incidence,  the  one-way  distance  traveled  by  compressional 
waves  in  any  layer  of  a K layer  media  is  the  thickness  of  the  layer.  For 
non-normal  incidence,  when  the  sensor  is  located  at  a distance  d from 
the  source  point,  the  one-way  distance  traveled  by  the  compressional 
waves  in  each  layer  is  a very  complicated  function  of  thicknesses, 
densities  and  velocities  of  aU  layers  (Ref.  5,  Figure  96,  for  example). 

The  total  distance  traveled  by  these  waves  in  all  K layers  is  a very 


i 

i 
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complicated  function  of  all  layer  thicknesses  and  all  angles  made  by 
each  section  of  the  wave  path  (Figure  22). 

Our  approach  is  to  approximate  the  one-way  primary  reflection 
paths  by  straight  lines  drawn  from  the  source  point  to  the  point  d/2 
at  the  bottom  of  each  layer  (Figure  23).  From  the  geometry  in 
Figure  23, 


T.  =1./  V. 
1 11 


( 1 ) 


so  that 


T.  = h.  / V. 
l li 


T.  = T.  ( l.  / h.  ) 
1 111 


( 2 ) 


( 3 ) 


l.  / h.  = 1 / cos  0.  = [1  +tan1 20.]^  = f 1 + (d/2)2/  £ h.  1^ 
ii  i i j=i  J 


( 4 ) 


which  can  also  be  written,  from  Eq.  (2),  as 


i.  / h.  = [ 1 + (d/2)2/  £ t.  V.  ]* 

1 i j=i  j y 


( 5 ) 


Our  final  expression  for  is  obtained  from  Eqs.  (3)  and  (5),  as 


d 

T.  = T. 

1 1 


. f 1 + d2/  4 £ t.  V.  J* 
1 j=l  J J 


( 6 ) 


Observe  that  t.  is  an  increasing  function  of  d,  that  for  d = 0, 
d i 

T.  = T.  , and,  that  becomes  closer  and  closer  to  t.  for  deep  layers, 
ii  l , l 

a 

We  propose  to  make  use  of  T.  in  the  following  manner.  Replace 

1 d 

Tj  in  all  of  our  earlier  models  by  t\  in  order  to  compute  the  signal 
observed  by  a sensor  which  is  offset  d units  from  the  source.  The 
implications  of  this  replacement  for  multiple  reflections  remains  to  be 
studied. 
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Figure  22.  Some  non-normal  incident  primary  reflection  paths 
taken  by  waves  in  reaching  a sensor  offset  d units 
from  the  source. 
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Figure  23.  Some  straight  line  approximations  to  non 
normal  incident  primary  reflection  paths. 
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VIII.  CONCLUSIONS 

In  this  paper  we  have  developed  state  space  models  for  a K layer 
media  system,  and,  have  shown  how  to  obtain  transfer  functions  from 
these  models.  Our  models  have  been  derived  by  examining  wave  effects 
within  a single  layer,  are  time-domain  models,  and  do  not  assume  equal 
travel  times  in  each  layer. 

Our  complete  model  is  of  order  2K;  it  generates  primaries  and  all 
multiples,  and  has  an  enormous  amount  of  internal  feedback  structure. 

Our  primaries  and  secondaries  models  are  each  of  order  2K;  the  former 
generates  just  the  primary  reflection  components  of  a seismogram,  where- 
as the  latter  generates  just  the  secondary  reflection  components  of  a 
seismogram.  Both  the  primaries  and  secondaries  models  have  no  internal 
feedback  structure. 

Our  development  of  the  primaries  and  secondaries  models  has  led 
us  to  a complete  (canonical)  decomposition  of  a seismogram  (see  Figure 
21).  The  utility  of  this  decomposition  remains  to  be  explored. 

Finally,  we  have  shown  how  our  normal  incidence  models  can  be 
extended  to  the  non-normal  incidence  situation. 

There  are  many  interesting  directions  for  future  work,  including: 

(1)  computational  aspects  of  the  various  models;  (2)  implications  and 

: 

applications  of  the  complete  (canonical)  decomposition  of  a seismogram; 

(3)  understanding  the  system  theoretic  aspects  of  causal  functional 
equations;  (4)  developing  optimal  state  estimators  for  our  models;  (5) 
studying  parameter  identification  problems  for  these  models,  etc. 

In  conclusion,  we  also  point  out  that,  whereas  we  have  presented 
all  of  our  results  in  the  context  of  a geophysical  model,  most  of  them 
(with  the  possible  exclusion  of  the  material  in  Section  VII)  are  applicable 
to  any  application  which  is  adequately  modeled  by  a lossless  wave 
equation  in  which  ray  theory  solutions  are  utilized.  We  shall  explore 
the  more  general  nature  of  our  results  in  a later  paper. 
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Figure  26,  First  residual  (i.e.f  complete  response  minus  the  primaries) 
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Figure  29.  Complete  response  sensed  at  normal  incidence  (lower  curve) 
500  ft.  (middle  curve),  and  1000  ft.  from  source.  The  top 
two  curves  were  obtained  by  means  of  our  approximate  non- 
normal incidence  transformations,  described  in  Section  VII. 


